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Abstract 

We introduce a mathematically rigorous analysis of a generalized spin- 
boson system for the treatment of a donor-acceptor (reactant-product) 
quantum system coupled to a thermal quantum noise. The donor /acceptor 
probability dynamics describes transport reactions in chemical processes 
in presence of a noisy environment - such as the electron transfer in a pho- 
tosynthetic reaction center. Besides being rigorous, our analysis has the 
advantages over previous ones that (1) we include a general, non energy- 
conserving system-environment interaction, and that (2) we allow for the 
donor or acceptor to consist of multiple energy levels lying closely together. 
We establish explicit expressions for the rates and the efficiency (final 
donor- acceptor population difference) of the reaction. In particular, we 
show that the rate increases for a multi- level acceptor, but the efficiency 
does not. 
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1 Introduction 



1.1 Transfer reactions and spin-boson model 

An important problem in chemistry and biology is to find electron transfer rates 
and transfer efficiencies in chemical reactions. A prominent example is the elec- 
tron transfer in proteins carrying out photosynthesis ([2-5]). The simplest re- 
actions are described by two states, a reactant (electron donor) and a product 
(electron acceptor). Before the reaction, the system is localized mainly in the 
reactant state, and after mainly in the product state. The passage from reactant 
to product is induced by two effects: a direct tunneling (hopping) and an indirect 
transition. The former originates from electron tunneling between reactant and 
product moieties, while the latter is due to the presence of thermal noise cre- 
ated by the many protein atoms and molecules in which the electron donor are 
acceptor are embedded. 

Denoting the reactant and product states by |R) and |P), respectively, a 
"Marcus model" Hamiltonian for the electron exchange has been used in [6,7], 



where Er and Ep are the reactant and product energies, and V is the direct 
tunneling constant. Both and Ep represent the collective energies of many 
particles (atoms and the molecules), corresponding to the reactant and product 
states of the protein environment. In the Marcus theory, the energy curves are 
taken to be harmonic in the collective position coordinate q, 



Here, / denotes the common force constant of reactant and product, gp is the 
equilibrium position of the product collective position (the reactant one being 
centered at the origin), and eo is the product-reactant energy difference. When 
describing the collective degrees of freedom of the environment quantum me- 
chanically, the reactant and product energies become the operators with the 
Hamiltonians of a collection of harmonic oscillators [7] 



The number of oscillators corresponds to the number of atoms in the proteins 
and is very large (of the order of 10 4 and more for the photosynthetic reaction 
center). Using these expression in the Hamiltonian i^Marcus above, we may write 



^Marcus = \R)E R {R\ + \P)E P (P\ + \R)V (P\ + |P)^(R|, 



E K = \fq\ E p = \f{q-qpf-e Q . 




H R V 
V Hp 
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which, under proper identification of parameters (see [7]), has the form of the 
spin-boson Hamiltonian (plus a constant term which we drop) 



is the bosonic "field operator" , expressed in terms of the creation and annihilation 
operators satisfying apo) a — a) a ap = if a ^ (3 and a a o) a — o) a a a = 1. 

One achievement of Xu and Schulten's work [7] is the identification of the 
Marcus theory model with a spin-boson system. This identification allows to 
take over results obtained previously for the spin-boson system. In particular, 
Xu and Schulten use Leggett et al.'s [8] expression for the transfer rate in the spin- 
boson system, and show that at high temperatures, it coincides with the transfer 
rate predicted by the Marcus theory. In this setting, it is assumed that the direct 
coupling V is very small. The term Va x is then viewed as a perturbation of the 
other terms in the spin-boson Hamiltonian (1.1). 

We recognize that the works cited above are of great importance in the field. 
Our aim for the present paper is to try to improve some points. 

1. Shortcoming. The derivation of the transfer rate for small V given in [8], 
and then used by [7], has two weak points. (1) The dominant term of 
the transfer rate is determined only heuristically (c.f. (3.31) of [8], errors 
in the perturbation expansion, i.e., terms of order V 2 and higher, cannot 
be estimated). (2) Additional approximations (of Born-Markov type) are 
made in the derivation. Their validity has not been (cannot be) verified 
rigorously (see before (3.32) of [8]). 

Remedy. We use the rigorous 'dynamical resonance method' [10,11] to 
find the dynamics of the reduced spin density matrix at all times t > 0, 
for arbitrary tunneling matrix elements V and arbitrary energy separations 
e. By rigorous, we mean that our results are derived by a mathematical 
perturbation theory in the coupling A between spin and bath, in which the 
remainder terms are controlled and are small uniformly for all times t > 0. 

2. Shortcoming. In the model (1.1) used by [7], the interaction term \a z ip(h) 
commutes with the main term |ea 2 (as V is assumed to be small). This 
means that in absence of direct hopping (V = 0), there is no electron 
transport at all. This is so since the populations are constant in time 
(diagonals of the density matrix in the energy representation). However, in 
reality, one would still expect electron transport due to the thermal noise, 
even if there is no direct hopping [12]. 



Hsb = Va x + \ea z + H R + \a z (p(h). 



(1.1) 



Here, a x and a z are Pauli matrices, and 




a 
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Remedy. We modify the Hamiltonian (1.1) by adding to the right side 
the term Xaa x ip(h) which induces population dynamics (electron transport) 
even if V = 0. Our analysis remains rigorous in presence of this term. 

3. Shortcoming. In the above model, both the electron donor and acceptor 
are assumed to have single levels. However, due to the complexity of the 
biological system at hand, it is reasonable to consider that either (or both) 
of donor and acceptor consist of a number of levels iVo, N\, centered around 
an average energy E D , E A . 

Remedy. We generalize the single-level system to the multi-level situation. 
We show that the multi-level model reduces to a single-level model with 

rescaled Hamiltonian matrix elements. 

While our approach allows for the above-mentioned improvements, we can 
only treat small values of the coupling between the donor-acceptor system and 
the thermal reservoir. Throughout this work, we assume that 



In this inequality, we consider A to be renormalized as to have the dimensionality 
of energy (i.e., (1.2) means that C|A| < — E&, where C is a small constant 
depending on the vairous parameters of the model, and which is such that the 
left side has the dimensionality of energy). 

The derivation of the electron transfer rate for small V in [8] is based on the 
fact that the spin-boson Hamiltonian, for V = 0, can be explicitly diagonalized 
(as in their model, the bath interaction commutes with the system Hamiltonian). 
As a result, the (heuristic) perturbation theory in V of [8] yields an expression for 
the electron transfer rate which contains all orders in the spin-reservoir coupling 
A. This approach cannot be carried out as soon as the interaction between spins 
and bosons is not proportional to o~ z , as the resulting Hamiltonian is not explicitly 
diagonalizable, even for V = 0. Instead of containing all orders in A and second 
order in V only, in our approach we obtain transfer rates (and the dynamics in 
general) to second order in A but to all orders in V, even in presence of a (de- 
) coherence altering interaction. We show in Section 5 that the order A 2 term of 
Leggett et al.'s transfer rate coincides with that obtained by our method. 

1.2 Transfer rates, separation 

Assume that initially, the donor is populated with probability one. In the course 
of time, the acceptor gains some population probability and for large times, the 
whole system converges to an asymptotic state. We call the transfer rate the speed 
at which the acceptor is populated, and the transfer separation the difference of 
acceptor minus donor population probability in the asymptotic state (after a 
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long time). The results presented here are immediate consequences of a much 
stronger result, Theorem 2.1 of Section 2, which gives the dynamics of the entire 
donor-acceptor density matrix, for all times. 



1.2.1 Single-level donor and acceptor 



The donor-acceptor-environment Hamiltonian is given by 



H = 



E D V 
V E A 



+ H R + X 



g D a 
a g A 



®<p(h). 



(1.3) 



Here, the first matrix is called the system Hamiltonian, it is the isolated donor- 
acceptor Hamiltonian, determined by the donor and acceptor energies, Ed and 
E A (with Ed > E A ), and the tunnelling matrix element, V G R. is the 
Hamiltonian of the uncoupled reservoir, a field of harmonic oscillators 



where we put H — 1, and a a , a) a are bosonic annihilation and creation operators. 
We take the oscillators to be in thermal equilibrium at temperature 1//3 > 0; 
they form a heat bath. An infinite- volume, or continuous-mode limit is taken in 
which the parameter a becomes the continuous boson momentum k G R 3 , see 
Section 4 for more detail. The coupling between the donor-acceptor system and 
the bosonic reservoir is described by the third part on the right side of (1.3). 
A G R is a coupling constant, <?d,<?a G R and a G R are interaction parameters 
responsible for energy conserving and energy exchange interactions. For a = 
(and when V = 0) we have the purely energy- conserving interaction (the situation 
considered in [7]). 

Key dynamical properties depend on a few system (donor-acceptor) and reser- 
voir (heat bath) quantities which we introduce now. The difference between the 
two eigenvalues of the system Hamiltonian is 




a 



n = V(e d - E A y + W 2 > o. 



(1.4) 



Define the parameter a, which can be positive, negative or zero, by 



V 



(1.5) 



a = 



Ed — E A 



The reservoir spectral density J(u>) is given, for uj > 0, by 



J(u) = y/n/2taiih(Pu/2) d(tu) + d(-uj) , 
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where C (uu) is the Fourier transform of 

C(t) = (e itH - V (h)e- itH ^(h))^ 



(1.6) 



the correlation function in the reservoir state at inverse temperature (3. See 
Section 4 for details. 

Transfer rates. The decay of the donor and acceptor populations, pd and 
Pa, is a complicated function of time in general (see Theorem 2.1). However, in 
the regimes where either Ed — Ea or V is small relative to the other one, the 
population decay (growth) is exponential in time, and we can identify a transfer 
rate. 



1. For \V\ « E D — Ea, we have 

PA(t) 



1 



1 + e~P Q 



+ 0(\ 2 + V), 



;i-7) 



where e is a complex (resonance) energy (see (2.1)). The remainder is 
uniform (homogeneous) in t > 0. The acceptor is populated monotonically 
exponentially, at the rate 7 re i ax = Imecb which has the form 



7rclax — 2 A 2 



a 



V 



9d - 9a 



Ed — Ea 

2i x4 



4aV 2 



coth(/3n/2)J(fi) 



+0(X V + A ). 
This result is obtained by taking a — > in Theorem 2.1. 
2. For E D - E A « |V|, we have 

PA(t) = \{l - Re e ite «) + 0(A 2 + E D - E A ), 



(1.8) 



(1.9) 



where is a complex (resonance) energy (see (2.2)). The remainder is 
uniform (homogeneous) in t > 0. The acceptor is populated exponentially, 
but modulated by cos(tRe £q), with rate 7r C i ax = Im£n, which has the form 



Orelax 



= AV^ 



(5-a + - 2a) 2 - (Ed - E A )(2 + ^—^-)(g A + 9d ~ 2a) 



+ (E D -E A ) 2 [ 



(9a - 9d + a) (9a + 9d - 2a) 



A 2 

T 



9a ~ 9d 



2V 2 

a(E D - E A ) 



V 
+ (1 + 



9a- 9d 
2V 



(Ed-Ea) 2 (9a-9dY 



2V J 4V 2 
x coth(f3Q/2)J(Q) + 0(X 2 (E D - E A ) 3 + A 4 ). 
This result is obtained by taking a — > oo in Theorem 2.1. 



2 ] C(0) 



(1.10) 
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Remark. The relaxation rates (1.16), (1.18) contain the product 

coth(/?ft/2) = \c(n) + c(-n) 

which depends on the temperature 1/(3 only via the bath correlation function C 
(see after (1.5)). We present the relaxation rates above involving the spectral 
density as this is customary in the literature. 

Transfer separation. We define the separation by 

S = pa(oo) -pd(oo). 

It measures how much the populations are separated after relaxation and satisfies 
— 1 < S < 1. The extreme cases S = ±1 correspond to complete localization of 
the final state in level one or two. S = means complete derealization (both 
levels equally probable). We obtain from Theorem 2.1 

1 2a 2 



S = -l + 



VI + 4a 2 1 1 + e-f™ 1 + VI + 4a 2 



(1.11) 



The separation does not depend on the initial state of the system. For Ed — Ex< 
< V (i.e., a large) we have S ~ 0, independently of the temperature. For 
V « Ed — E^ (a small) we have S w — 1 + 1+c 2 -gsi = tanh(/3f2/2), which becomes 
S at high temperatures (f3 — > 0), and S ps 1 at low temperatures ((5 — > oo). 
We conclude: 

For large hopping constant (Ed — E\ « V) the donor and acceptor are pop- 
ulated equally in the long run (S ~ 0). The same happens for small hopping 
constant (V <K Ed — Ex) at high temperature. However, for small hopping con- 
stant (V « Ed — Ex) and low temperature, the acceptor is fully populated in the 
long run (and the donor has probability zero). 

1.2.2 Multi-level acceptor model 

For an A^-fold degenerate acceptor, the total Hamiltonian is 





E D 


V ■■ 


■ V 




9b 


a ■ 


a 




V 


E A 




+ H R + X 


a 


9A 




H = 












V 




E A 




a 




9A 



®<p(h). (1.12) 



We introduce the donor state, <^d, and the collective acceptor state, ax, 



1 





OA 





1 



(1.13) 
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In the basis {<^d,o"a}, (1-12) takes the form 



H = 



V^Nl E A 



+ H K + \ 



ay/W 9a 



<p(h). 



(1.14) 



This Hamiltonian is of the form (1.3), but in a different basis, and with rescaled 
off-diagonal matrix elements. We can thus take over results obtained for the two- 
level model (see Section 2 for details). As explained below, the final populations 
depend on the initial density matrix in the multi-level case. In this section, we 
assume that initially, the donor is fully populated and there is no donor-acceptor 
entanglement initially. In other words, the initial donor-acceptor density matrix 
is \(pv)((pn\- 

Remark. We consider here V independent of N A , so the direct do nor- accep- 
tor interaction is of the size VN A . One may compensate this growth by scaling 
V with a negative power of N A . We do not pursue this question in the present 
manuscript. 

Transfer rates, multi-level acceptor model. As in the case of the two- 
level system, the decay of populations is exponential in two limiting cases. 



1. For \/Na\V\ « E D - E A , we have 

1 1 " 

PA(t) 



3 iteo 



N A 1 + e~P n 



+ 0(\ 2 + V 2 ), 



1.15) 



where e$ is a complex resonance energy, depending on N A (see Theorem 2.3). 
The acceptor is thus populated exponentially quickly and monotonically, at 
the rate 



Orclax 



2N A X 2 



V 



gp ~ 9a 
Ed — E A 



AN A a 2 V 2 



+0(X 2 V 3 + X 4 ). 
A proof of this is obtained by taking a small in Theorem 2.2. 
2. For E D - E A « y/W\V\, we have 



coth(/5f2/2)J(fi) 

(1.16) 



PA(t) 



2N A 



■(1 — Re e iten ) + 0(\ 2 + E-q — E A ), 



(1.17) 



where Eq is a complex resonance energy, depending on (see Theorem 
2.3). The acceptor is populated exponentially, modulated by cos(tRe Eq), 
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with rate 

%elax = 



A 



(9a + 9b ~ 2a\fN A ) 



-(E D - E A )(2 + ^—^r)(g A + 9D - 2a y/K) 



+ (E D -E A ) 2 [ 



V^N A 

(9a - 9d + ay/N A )(g A + g D - 2ay/W) 



2V 2 N A 



+ (1 + 



9A-gv , 2 
2V^N A ~ 



f] C(0) 



+- 



9a- 9b- 



a(E D -E A )\ 2 (E D - E A ) 2 (g A - g D ) 



2V 



4V 2 N A 



x 



coth(/3ft/2)J(ft) + 0(\ 2 (E B - E A ) 3 + A 4 ). 



;i.i8) 



A proof of this is obtained by taking a small in Theorem 2.2. For large N A , 
we obtain 



7 ; elax ~ 2V27r\ 2 a 2 N A C(0) 



9a - gD - 



a(E D - E A ] 



2V 



(Note that Q « 2V\fN A in this case.) 



J{2V^N A ). 

(1.19) 



Transfer separation for the multi-level system. Since all acceptor levels 
are populated equally at all times, we define the transfer separation for the multi- 
level acceptor system by 

S = N aPa (oc)- Pd (oc). (1.20) 
Then S = 2 (\cr A ) (a A \) ^ — 1, which is given again by (1.11), but with 

a = ^N~ A V . (1.21) 

For large N A , we have 5 = 0, which means that all donor levels share probability 
1/2 and all acceptor levels share probability 1/2 as well. 

Discussion: effects of multiple levels. (1) Transfer rate and degener- 
acy. The effect of the A^-fold acceptor degereracy is to multiply the hopping 
coefficient V and the decoherence coefficient a by \/N A , while it does not affect 
energy conserving parameters (see (1.14)). Accordingly, the relaxation rates are 
(roughly) acquiring a factor N A (as they are proportional to the square of the co- 
efficients V and a), see (1.16) and (1.19). Therefore, N A -fold acceptor degeneracy 
speeds up the transfer process, the rate being proportional to N A . 
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(2) Asymptotic population and degeneracy. Due to the scaling, the Hamilto- 
nian (1.14) becomes, for large N A , 



VJN A 



1 

1 



1 

1 



®(p(h). 



Therefore, the off-diagonal density matrix elements of the donor-acceptor system 
(in the basis {</?d,ca}) is time- independent. Moreover, the system approaches 
the Gibbs equilibrium state in the long run, and the latter is, up to 0(A 2 )-terms, 
equal to 



This is why, in presence of many acceptor levels, the final populations both in the 
donor, and in all acceptors together, equal one half each. Hence pa = 1/(2 A^). 

(3) Separation and degeneracy. In line with equal total distribution of donor 
and acceptor levels (previous point), the separation must vanish for large Na- 

(4) Degenerate donor. If the donor is Afo-fold degenerate and the acceptor 
is simple, then all the above formulas for the transfer rates and separations are 
the same, upon replacing Na by A^d- 

(5) Both acceptor and donor degenerate. Some of our results hold if both 
the donor and the acceptor are degenerate. However, so far, we have not been 
able to find the dynamics of both the donor and acceptor in this setting, see the 
explanations in Section 2. However, a consideration as in point (2) above gives 
the following asymptotic result: if the donor and acceptor have degeneracies N D 
and Na, respectively, then the transfer rates scale as ^N^Na, and the separa- 
tion becomes zero for large Nd and Na- Asymptotically, each donor level has 
probability 1/(2A^d) and each acceptor level has probability 1/(2A^). 

(6) Quasi-degenerate levels. Our approach is also applicable if the levels are 
not exactly degenerate, but, say, spread around average values E-q and Ea- More 
precisely, if, similar and in addition to (1.2), the spread AE of the levels satisfies 
AE E D — Ea, then the formulas below for transfer rates and separation give 
the correct lowest order terms in AE. In principle, one can calculate corrections 
of order AE, but this is rather complicated. 

(7) Dependence on initial condition. Consider a doubly-degenerate acceptor, 
system (1.12) with 3x3 matrices. This system has two invariant states (the kernel 
of H has dimension two). One is immediately seen to be 



1 



T = 



y/2 




1 



^1, 



where is the equilibrium state of the reservoir (satisfying HrQh = 0). The 
other stationary state is given by the 2D-reduced Gibbs equilibrium state of 



10 



the two-dimensional system interacting with the reservoir, expressed in the ba- 
sis {(/?d ? 0"a}) (1-14). After tracing out the reservoir degrees of freedom, and 
to lowest order in the system-reservoir interaction A, this 2D-Gibbs state is 
f2s oc exp(— /3H S ), where H s is the first matrix on the right side of (1.12). For 
instance, if V = 0, then 

O e-^|v9 D )(v9 D |+e-^|a- A )(a A | ( , 
" s = e -^D +e -^A • ( h22 > 

Note that in the original basis {y?D, fi, ¥2}, in which the matrices in (1.12) are 
presented, the final state (1.22) is not even diagonal. 

A general initial condition will converge, for large times, to a superposition of 
the two invariant states. If the initial condition belongs to one of the invariant 
subspaces span{y9D,aA} or Cr, then its final state will be the 2D-reduced Gibbs 
state or r. In the above results, we assume that the initial condition is the pure 
state entirely concentrated on the donor. It belongs to the subspace spanned 
by {(/?d,oa} and hence the final state is the associated Gibbs state. A similar 
dependence of the asymptotic state on the initial condition holds for any acceptor 
dimension. We note that a deviation of the Gibbs equilibrium as a final state of 
an open system has also been observed in [14], in the setting of a spin coupled 
symmetrically to a bath of other spins in thermal equilibrium ('spin star system'). 

(8) Continuous acceptor, or sink, or Wigner-Weiskopf limit in presence of a 
heat bath? In so-called "sink" or Wigner-Weiskopf models, a dissipative part of 
the system is modeled by considering a system of interest (the donor) coupled to 
N energy levels (the acceptors), and the limit N — > 00 is taken in order to obtain 
irreversible phenomena (decay). See, for example [13], p. 36 and following. In 
taking the continuum limit, the typical spacing between individual levels, AE, 
is decreased more and more and at the end the levels are characterized by a 
continuous density of states. The typical role of a sink is to depopulate the 
donor exponentially quickly (at a rate proportional to the density of states). One 
may ask if, in taking N\ — > 00 in our model, we obtain the same dynamical or 
asymptotic results as for a sink model. The answer is negative for the reason we 
explain below. However, we point out that the rigorous derivation of the Wigner- 
Weiskopf model has never been done (to our knowledge) in the presence of an 
additional heat bath. 

The following transition occurs when taking the AE — > limit in presence 
of a thermal bath. Let A be the coupling strength of the donor-acceptor sys- 
tem with the thermal bath. If A AE, then the donor-acceptor energy levels 
"are well defined" (as without interaction with the heat bath). In this case, 
the asymptotic state is given by the Gibbs state of the donor-acceptor system 
(oc exp(— j3 Hs), modulo corrections small in A). The final donor population is 
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consequently (consider V = in (1.12)) 



Pb = 



e -l3E u _|_ ]\[ e -l3E A ' 



which decreases as 1/N. In this regime the transfer has a good efficiency (almost 
total depopulation of the donor), as is the case in the Wigner-Weiskopf model. 
However, as we increase N, we decrease AE and we reach the regime A ys> AE. 
This case is a perturbation of the totally degenerate situation (AE = 0) which 
we treat in the present paper. The 2D- reduction (1.14) then takes place, and we 
obtain a final donor population (see after (1.21)) 



which means that the transfer is not very efficient (and cannot be made more so 
by increasing the number N of acceptors). 

To sum up: Our analysis holds for arbitrary N, but in the limit N — > oo, 
it does not give depopulation of the donor. The depopulation however holds in 
sink-models without a thermal bath. 

2 Main results: details 

2.0.3 Two-level donor-acceptor model 

The total Hamiltonian is given by (1.3). We point out that J(uj) is independent 
of the temperature. It has the explicit representation (u > 0) 



where the integral is over the two-dimensional sphere, and the function h = h(k) 
is written in spherical coordinates for k = (oo, E) e M 3 . (In the last expression 
for J, the oscillators are indexed by a — > k in the continuum, or infinite volume 
limit for the momentum k G 1R 3 .) 

Let p t be the reduced donor- acceptor density matrix, when the degrees of 
freedom of the environment are traced over. We assume that initially, the entire 
system is in a state of the form 



where p is an arbitrary two-state initial density matrix, and pn is the initial 
state of the reservoir, taken to be at equilibrium at inverse temperature (3. Then 
we have the dynamical equation 



Pd ~ 1/2, 




Pin = Po ® PR, 



Pt = Tr R [e- itH p- m j tH ] 
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for the reduced density matrix. Here, the trace is taken over the reservoir space. 
The two eigenvalues of the donor-acceptor Hamiltonian are 

X - [E B + E a ± v/(£ d - E A ) 2 + 4V 2 , 

resulting in the eigenvalue difference (1.4). We define the 'resonance energies' by 

e = 2iA 2 [a + ( ^ A ~^ D)a]2 coth(^/2)J^) (2.1) 
1 + 4cr 

£n = e /2-n-X 2 X + i\ 2 y^/2 Y 2 C(0) (2.2) 

where 

4a 2 4a 

Y = <* - ^ VTTWWTTW^) + 2JD - a VTTW + * - Ed - (2 ' 3) 

and 

(8ck \ f°° 

4g D - a ^——^ - Yj 1m jf (p(t)^ df 

+2^ [a+ [ 9A + ~^ )a]2 Re sin(fii) (^(i)^ df . (2.4) 

Here, ( ) /3 denotes the average in the thermal equilibrium of the environment. We 
denote by [pt]ij the matrix elements of the reduced density matrix of the two-level 
system in the basis {[1 0]*, [0 1]*}. The following is our main result. It describes 
the population dynamics of the two-level system, identifying a main part and a 
remainder, which is of order 0(A 2 ), homogeneously in (independent of) time. 

Theorem 2.1 The dynamics of the product (acceptor) probability [pt\22, for t > 
0, is given by 

N22 = N 2 2 TT ^{e ite0 + 4a 2 Ree ite "} 

-Re([p ]i2) TT ^{e lte »-Ree^} 
2a 

+Im(Mi2) IT ^Im e 1! 



,iten 



[1 + e-^]" 1 ^ _ e ite ° + VTT4c7 2a 2 



VI + 4a 2 v ' 1 + VI + 4a 2 I + 4a 2 

2a 2 

~ -Re e iten +0(A 2 ). (2.5) 



1 + 4a 2 



The remainder term 0(A 2 ) is independent of t > (and of E D , E^,V,a varying 
in bounded sets). 



13 



Remark. In the "usual" setup [9-11] for the derivation of the reduced density 
matrix (2.5), we start with a diagonal system Hamiltonian. Then the diagonal 
density matrix elements evolve jointly, and only the 'resonances bifurcating out 
of the zero eigenvalue' (eo here) are present in their evolution. However, in the 
present setup, H$ is not diagonal, and as a result, the evolution of the diagonal 
involves the initial condition of the off-diagonal density matrix elements, and 
the evolution also depends on the resonances bifurcating out of the non-zero 
eigenvalues {eq here). 



2.0.4 Multi-level donor-acceptor model 



We consider an iV D -fold donor and an A^.-fold acceptor with energies and 
E\, respectively. The energy levels may be distributed around these two fixed 
energies, provided their spread is small. The total Hamiltonian is 



H = 



Ejy 




V . 


. V " 




E B 


V . 


. V 


V . 


. V 


E A 




V . 


. V 







+A 



' 9d 




a . 


a 




9d 


a . 


a 


a 


a 


9a 




a 


a 




9A . 



® (fi(h). 



(2.6) 



The donor- acceptor space is partitioned into iVo levels Ed and Na levels Ea- By 
introducing the vectors 



0"D 



1 





the Hamiltonian (2.6) can be written as 



0~A 






1 



(2.7) 



H 



£ + VViV D iV A {|(x D >(o- A | + |oa)(<7 D |} + # r 

+A|^ + a v / iVDA^(k D )(oA| + |<7 A )(<7 D |)}®</?(/l), 



(2.8) 
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where 



E = dia,g(E D ,...,E D ,E A ,...,E A ), 
Q = diag(0D,...,0D,£k,...,flO. 



(2.9) 
(2.10) 



Reduction to two-level system. From (2.8)-(2.10) we see that H leaves 
the subspace 

H = span{o- D , a A }® Hr 

invariant. This means that if x G 'H then H\ £ It is the same as saying that 
H can be written as a block-diagonal matrix in the decomposition 

H = H@H ± . 

Of course, each block is still infinite-dimensional due to the reservoir degrees of 
freedom, but the block in H involves only two donor-acceptor vectors, namely 
o"d and a A . This is why the the multi-level donor-acceptor model has a two- level 
formulation. In the basis {(Td,o"a}, H in (2.6) takes the form 



H 



E D _ Vy/l^W 

E A 



+H R +\ 



9d ay/N D N A 
aVN D N A g A 



)<p(h). (2.11) 



This Hamiltonian is of the form (1.3) with rescaled off-diagonal coefficients V — > 
VVN D N A and a ->■ ay/N B N A . 

Symmetry. Due to the symmetry of the Hamiltonian H, (2.8), the donor- 
acceptor density matrix has a special structure. For 1 < i, j < N D + N A , let Uij be 
the unitary operator wich exchanges labels % and j. In other words, Uijipi = tpj, 
Uijipj = (fi and Uijipk — (fik if k ^ i,j, where {(pk}^=i NA is the energy basis (in 
which (2.6) and (2.7) are expressed). 

We consider initial density matrices p which are symmetric with respect to 
permutation within the donor and within the acceptor subspaces, 



Po 



(2.12) 



if 1 < z , j < iV D and if iV D + 1 < i, j ' < N D + N A . An example of a symmetric 
initial state is 

p = -^diag(l,...,l,0,...,0), (2.13) 

iVD 

in which each donor degree of freedom is populated equally likely. 

Theorem 2.2 (Symmetry) Suppose that the initial state is symmetric as in 
(2.12). The reduced density matrix of the donor- acceptor system has the form 



Pt = 





* 


* 


x A (t) 
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where the * represent some matrices, and where Xn and Xa are square matrices 
of size N D and Na, of the form 



X D (t) 



Pd %d 
xd Pd 

x D ■■■ 



xb 



■ ■ x D 

xd Pd 



x A (t) 



Pa x a 
xa Pa 

xa ■■■ 



xa 



■ ■ xa 

xa Pa 



The off-diagonal matrix elements of Xd are all equal and real, and so are those 
of Xa- All diagonals of X D are equal, and so are those of Xa, and they satisfy 
PaN a +PbN d = 1. 



Proof of Theorem 2.2. For 1 < i,j,k < N D we have 

\pt\ij = Ti{U kj U kjP ^ tH \ Vj )(ip t \e- itH ) 

= Ti( P ^ tH U k] \ V] ){ Vl \U k ^ itH ). (2.14) 

If % = j then U k j\(pj) (tpi\U k j = \^ k ){^ k \ and (2.14) means that [p t }n = [p t ]kk- If 
% ^ j then: U k jtfi = tfj if k = i (so (2.14) gives [p t ]ij = \pt\ji) and U kj <pi = ipi 
if k 7^ % (in which case (2.14) gives [pt]ij = \fit\ik)- This shows that A is of 
the form as given in the theorem. Repeating the same argument for indices 
iV D + 1 < i, j, k < N D + Na yields the form of B. The relation between p and q 
is obtained from Tr p t — 1. ■ 

Consider the average (\(Ta) {<ja\) t , where we recall that <ta is given in (2.7). On 
the one hand, we have 

(\a A )(a A \) t = Tr(p e itH \a A )(o- A \e- itH ) 

1 N A +N D 

- W E Tr( Po e^|^)(^|e-^) 

1 N A +N D 

= Na E ^ 

i,j=N n +l 

PaN a + x a (NI-N a ) 
= t- t =pa + x a {N a - 1), (2.15) 

where pa and xa are the matrix elements of the acceptor block given in Theorem 
2.2. On the other hand, by comparing (2.11) with (1.3), (\o'a) (ca\) t equals \ftt\22 
in the formalism of the two-level model. Theorem 2.1 thus yields the following 
result. 
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Theorem 2.3 Suppose that initially, the donor degrees of freedom are populated 
only, with equal probability l/N D , as in (2.13). Then we have 

[i + e-^]-i e iteo + Vl + 4a 2 2a 2 

2 ' r r Re e ite » +0(A 2 ). (2.16) 



1 +4a 2 " 

i/ere, a, Q, e , Eq are given as in (1.4)-(2.4), but with a and V replaced by a^NvN A 
and Vy/N D N A , respectively. 

By proceeding in the same way, one finds 

(\o-D)(o-D\) t = PD + x D (N D -l), (2.17) 

where pn, xb define the donor block Xd defined in Theorem 2.2, and (\o~d) (o~d\) t 
is equal to \fit]ii = 1 — \ftt\22 in the two-site model (Theorem 2.1). This means 
that (\cr D )(a D \) t + (\<t a ) (<t a \) t = 1, hence (2.15) and (2.17) give 

Pa + x A (N A - 1) +p D + x D (N D - 1) = 1. (2.18) 

Together with the equation pd^d + PaN a = 1 (see Theorem 2.2), (2.15), (2.17) 
and (2.18) are four equations for four unknowns p A ,pn,x A ,Xr>- However, only 
three of those equations are independent (as (2.18) is the sum of (2.15) and 
(2.17)), so the solution is indetermined. In case N A = 1 or Nd = 1, one variable 
is eliminated (for instance, if Nd = 1 then x-q is not present), and the system of 
equations can be solved. We consider this next. 

Single-level donor to multi-level acceptor. We look at N D = 1 and 
A^ > 1 arbitrary. Then 

p D = 1 - (|<7 A )(<7 A |) t , p A = -^-(\o- A )(a A \) t , x A = ^-(\a A )(a A \) t . 

The first equation comes from (2.17) and (|<7d) (&i)\) t + (\cr A )(a A \) t = 1. The 
second equation then follows from pt>Nd +p A N A = 1. Finally, the third equation 
comes from (2.18). 

Multi-level donor to single-level acceptor. Here N D > 1 and N A = 1. 

Then 

_ 1 ~ (\0-A)(0-A\) t -/UW^IX T _ I" (|0-A)(0A|) t 

Pd 77 , Pa - {\o- A ){cr A \) t , x D . 

iVD iVD 

The statements about the transfer rates and separation in Secton (1.2.2) follow 
directly from the above formulas. 



17 



3 Proof of Theorem 2.1 



The dynamical resonance method [10,11] gives the reduced dynamics of the spin 
system (donor-acceptor system). Its starting point is a diagonal system Hamil- 
tonian. Then a rigorous perturbation theory is applied, tracing out the reservoir 
degrees of freedom and describing how the system energies become complex ('res- 
onance energies') and lead to decay. Those complex energies are the eigenvalues 
of a (non-hermitian) effective energy operator, called a "level shift operator" . The 
task is thus to diagonalize the 2x2 system Hamiltonian, and then to calculate 
and diagonalize the level shift operators. We do not carry out all details as this 
would take up too much space. 
The system Hamiltonian 



E D V 
V E A 



is diagonalized by the unitary 



U 



v 



_Ci 



V C2 



(3.1) 



where 



Ci,2 = \{Eb + e aT V(E D -E A y + w 2 }. 

In the new basis, the Hamiltonian (1.3) has the form 

H = UHU- 1 = H s + H R + XW <g> <p(h), 



where 



i/s:=^sf/" 1 = diag(Ci,C2) 



and W = UWU~ l . As the system Hamiltonian is now diagonal, we can apply the 
dynamical resonance theory [10, 11] to find the dynamics of the reduced system 
density matrix. We outline the most important steps. 

In the Gelfand-Naimark- Segal Hilbert space representation of the system, the 
density matrices on the spin-boson Hilbert space, C 2 ® J r (L 2 (M 3 , da; 3 ), are identi- 
fied with vectors of a new Hilbert space "Hqns = C 2 (g>C 2 (g)J : '(L 2 (]Rx S* 2 , du x dS). 
This "doubling" of the space is described in detail in [10,11]. Here, J-'(X) is the 
Fock space over the one-particle space X [1]. For I = L 2 (l x S 2 ,du x dS) it 
carries the creation operators and annihilation operators a(u, S), a){u, S) satis- 
fying [a(u, S), o)(u\ £')] = 5(u — «')£(£ — £') (Kronecker deltas). The Liouville 
operator is defined by 

K = L s + L R + XI, 
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where L$ = H$ <S> Us — Is ® i?s, £r = dr(u) is the second quantization of 
the operator of multiplication by the argument u G 1 in L 2 (I x S* 2 ). The 
operator / represents^the interaction between the spin and the bosons. The 
explicit form is / = W ® t s <g> <pp(h) - JA 1/2 (W <g> l s ® (f/3(h))JA 1/2 , where 
J, A are the modular conjugation and the modular operator associated to the 
vector fig ® Or € "Hqns- Here, £7$ is the trace state of the spin and Or is 
the vacuum vector, representing the equilibrium state at temperature 1/(3 > 
of the infinitely extended bose gas. The modular data is a concept of Tomita- 
Takesaki theory of von-Neumann algebras. We do not define these objects here, 
but refer to [1,10,11] for details. The spectrum of L s consists of energy differences 
of H s . The spectrum of L R has a simple eigenvalue zero with eigenvector R 
and is otherwise continuous, covering the whole real axis. The four eigenvalues 
of K = L$ + Lr. are ±(Ci — C2) (each simple) and zero (twice degenerate). 
These eigenvalues are embedded in the continuous spectrum. The achievement 
of the dynamical resonance theory is to describe the instability of the eigenvalues 
under the perturbation XI. More precisely, these eigenvalues become the complex 
eigenvalues of a spectrally deformed version of the operator K. To each of the 
eigenvalue e of K$ is associated a "level shift operator" A e . The eigenvalues of A e 
are (up to fourth-order corrections in A) the complex eigenvalues. For e = ±(Ci — 
(2), A e is simply a number, namely the operator acting on the one-dimensional 
space Cy?i ® y?2 <E> Or and Cif2 CS> y^i ® Or, where fi^ is the canonical basis in 
which Hs is diagonal. The operator A is two-dimensional (as zero is doubly 
degenerate). We now give the explicit form of these level shift operators. (Their 
definition and calculation is rather straightforward, albeit somewhat lengthy - 
we refer to [10, 11] for a general formulas for the level shift operators). We have 



where O is given in (1.4) and where C is the Fourier transform of the correlation 
function (1.6) (see also Section 4). This operator is written in the basis {ipi <S> 
</?i, (/9 2 <8y?2}- The eigenvalues of (3.2) are and e (see (2.1)). The eigenprojection 
onto the eigenvalue zero is 




(3.2) 




while that on the eigenvalue e is 



Q ( 



(i) 



c(o) + c(-o) L 



-1 



1 

-1 




-C(O) 

6(-n) 
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We point out the formula C(Q) + C(—Q) = fi 2 coth(/3Q/2), see Section 4. Simi- 
larly, one finds for the level shift operator associated to ( 2 — Ci — ^ 

An = £n \<f2 <E> <fi)(ip2 <S> (pi\, (3.3) 

where en is given in (2.2). Similarly, = — |y?i ® </? 2 ) (v^i ® V^l- 

This information is sufficient to give the dynamics of the reduced system 
density matrix p t in the basis in which H s is diagonal (see e.g. Theorem 2.1 
in [10] or Theorem 2.1 in [9]). For example, 

[Pt]u ■= (<Pi,Ps<Pi) = A(H; H)[po]ii + A(H; 22)[p ] 22 + 0(A 2 ), (3.4) 

where 

A t (ll; kk) = (ip k ® (p k , Qo°Vi ®<Pi)+ e lte ° ' (<fk ® (fk, Qo^i <8> • 

The remainder term in (3.4) is uniform in t > 0. To obtain the dynamics of the 
reduced system state in the original basis (in which H s is not diagonal), we undo 
the base-change implemented by U, according to p t = U^ptU, c.f. (3.1). This 
yields (after some algebra) the relation (2.5) for [p t ] 2 2 = 1 — [pt]n- • 



4 Relation between bath correlation and spec- 
tral density functions 

In this section, we derive the relation between the spectral density and the bath 
correlation functions, see (4.13). As the spectral density function is defined in 
the physics literature for environments of harmonic oscillators labelled by a dis- 
crete parameter (momentum), we first identify all quantities of our (continuous 
momentum, i.e., infinite volume) model with the equivalent discrete counterparts. 

In the seminal paper [8] on the spin-Boson system, the following model is 
considered. A spin is coupled to a bath of oscillators. The Hamiltonian of a single 
Hamiltonian, labelled by a, is H R , as introduced before (1.1). The interaction of 
the bath oscillators with the spin is given by 

^qoo- z ^2c a x a , (4.1) 

a 

where a z is the Pauli operator, 

- 1 / u 

%a — /jr ( a a + a a)i 

y/2m a UJ a 

see equation (1.4) and p. 7 of [8]. Here, q is a coupling constant, c a are real 
numbers and a a ,a^ are annihilation and creation operators, satisfying a a a^ — 



20 



a^a a = 5 a/ 3 (Kronecker symbol). Leggett et al. define the spectral density in [8], 
equation (1.5), by 

^> = fE^(— «>■ 

a 

Our interaction in (1.3) has the form 



(4.2) 



A 



9d a 
a 9a 



<p{h). 



(4.3) 



This interaction is identified with the quantities in Leggett et al.'s work by making 
a discretization of momentum space of the free field in a box of (very large) side 
length L: 



d 3 k 



¥) E 



2tt 



-3/2 



Then, taking the function h(k) in (4.3) to be real valued, we have 
i r 1 / 9tt \ 

m = 71 L m{a ^ k) + a(A;)}d3A; ~ 7f (tJ ^ (4 + «*)■ 



Thus our interaction (4.3) is of Leggett et al.'s form (4.1), under the following 
identifications: 



9d = ~9a 
A 

2vr\ 3/2 
L 



V2 





1/2 



E = E 

1 - 2?T r77Q ^ 



(4.4) 
(4.5) 
(4.6) 

(4.7) 
(4.8) 



We evaluate the Fourier transform of the correlation function, 

1 Z" 00 

C(u>) = -j=J e-"*(<p(t)<p)dt. 
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We have 

a,a' 

= ^^(e^V + e-^K + l)), (4.10) 

a 

where n a = (a^a a ) is the average occupation of mode a. Taking the Fourier 
transform of (4.10), and using that 

/oo 
e i{ul ~ UJ ' )t dt = 2n5(uj-u'), (4.11) 
-oo 

we obtain, for wGl, 

C(u>) + d(-ou) = v^72 ^(1 + 2n a )5(cu - u a ). (4.12) 

a 

For a thermal reservoir, we have n a = (e^ — l) -1 , so 

C(u) + d(-co) = y/2/^coth(pu/2)J(uj), (4.13) 
where J{oj) is the spectral density (4.2). 



5 Recovering Leggett et al.'s results for small 
tunneling 

In the setting of the Marcus theory of electron transport [6, 7] the system- 
envoronment interaction is diagonal, a = in (1.3), and the tunneling element 
V in (1.3) is considered to be very small. Transfer rates are then obtained per- 
turbatively in orders of V. In this section, we show that the rate obtained in our 
paper is the same expression as that obtained in Leggett et al.'s work [8]. 

For a = and \ V\ « E D — E A and to lowest order in A and V, our relaxation 
rate (1.8) is 

W = 2AV 2 ( g D ~ | ) coth ((3(E B - E A )/2)j(E D - E A ). (5.1) 

\ tin — tiA J 

Since Q = E& — E A + 0(V), we have replaced Q in (1.8) by the energy difference 
Ej)—E A in (5.1). To compare our rate with Leggett's, we first harmonize notation. 
Our parameters are identified with those of the spin-boson model in [8] according 
to (4.4)-(4.8) and 

V = -A/2 (5.2) 
E D -E A = e, (5.3) 
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where A and e are the energy gap and the tunnelling constant of [8]. Therefore, 
in Leggett's notation, our transfer rate (5.1) is given by 



7reiax = f ^coth(/3e/2)J(e). (5.4) 

Leggett et al. use the heuristic "golden rule" approach (Section III.D of 
, which is a formal perturbation theory in the tunneling constant A (plus 
some further approximations), to derive the following expression for the electron 
transfer rate (see equation (3.38) of [8]) 

r- x = A 2 / dtcos(et)cos[(g 2 /7r)g 1 (t)]exp-[(g 2 /7r)g 2 (t)], (5.5) 
Jo 

where 

QJt) = [ ^Msinwtdw, (5.6) 
Jo w 

^ / s f°° J(w)(l — cos cut) , . . , . „ s 

Q 2 (t) = / — = coth(/3u;/2)dw. (5.7) 

Jo u 

Note that (5.5) contains all orders in qo, but our result only yields the second order 
(i.e., the q 2 = A 2 term). While we use a mathematically rigorous perturbation 
theory, Leggett et al. proceed as follows. First, they apply a suitable unitary 
transformation to the spin-boson Hamiltonian, transforming it into ( [8], equation 
(3.30)) 

H' = -\A(a + e~ in + H.c.) + \ea z + £ ^(m a u 2 a x 2 a + p 2 Jm a ). (5.8) 



a 



Here, a + = \{cr x + icj/) and Q = J2 a {QoCa/ 'm>aWa)Pa- For A = this operator is 
explicitly diagonalized. Then they use perturbation theory in (5.8) for small A, 
and obtain (5.5) to second order in A. Note that in applying perturbation theory 
to (5.8), the perturbation is — |A(a + e -1 ^ + H.c), which depends on q , indirectly 
via Q. Even though this perturbation is bounded by const. |A|, independently of 
q , the approximations made to arrive at (5.5) are complicated, and it is not clear 
if (5.5) is the correct expression for relatively large go- (The correction term could 
be, for instance, of order A 4 g 2 , which would then dominate the main part, (5.5), 
for A 2 gp > 1-) For small values of q , the perturbation series can be controlled, 
and formula (5.5) should be rigorously correct. 
Expanding (5.5) for small q , we have 



■t = A'/ it e "' + e "" 



o 



i - ^W) 

71 



+ O(AY ). (5.9) 
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Using equation (5.7) we can integrate explicitly over the variable t (c.f. (4.11)), 



Comparing (5.10) and (5.4) shows that: 

The expression for the transfer rate obtained by the resonance method, to sec- 
ond order in the tunneling matrix element (V ) and second order in the interaction 
between the spin and the reservoir (\), is the same as that quantity obtained by 
the heuristic "golden rule" method used in Leggett et al. [8]. Note that the equality 
of the two expresssions holds for arbitrary spectral density functions J. 

Remarks. (1) [8] gives an expression for the transfer rate which contains all 
orders in the coupling q to the reservoir, see (5.5). However, its derivation is not 
controlled: it is not known if the 'remainders' in the perturbation arguments are 
smaller than the 'main terms'. Our method is mathematically rigorous, giving 
bounds on all remainder terms in the perturbation theory arguments, in the 
parameter regime where the coupling to the reservoir is small enough. This means 
that our main term is guaranteed to be larger than the perturbation corrections. 
Higher order terms in the spin-reservoir coupling can be calculated rigorously 
with the resonance method, but their derivation is lengthy, and we did not check 
if they coincide with Leggett et al.'s expressions. 

(2) It has been shown in [7] that the spin-boson system (the same as in [8]) 
gives the same reaction rate for large temperatures (in the physiological regime). 
Since we have just shown that our results coincide with those of [8], we have that 
the resonance theory produces the same transfer rate as the Marcus theory, in the 
mathematically controllable parameter regime (small tunnelling matrix element) 
and at physiological temperatures. 
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